library(RANN)
library(arrow)
library(dplyr)
library(tidyr)
library(cluster)
library(circlize)
library(ComplexHeatmap)
library(SingleCellExperiment)
source("code/_utils.R")ctx
preamble
loading
sce <- readRDS("outs/fil.rds")
ist <- readRDS("outs/lv1.rds")
jst <- readRDS("outs/lv2.rds")labeling
sce$lv1 <- (kid <- ist$clust)[match(colnames(sce), names(kid))]
kid <- unlist(unname(lapply(jst, \(.) .$clust)))
idx <- match(colnames(sce), names(kid))
table(sce$lv2 <- kid[idx], exclude=NULL)
BEC epi fib FRCcts FRCpv FRCtcz LEC mCAF AML DC
1034 212 653 941 744 3753 945 2016 2943 1700
mac.pi mono.c mono.nc pDC TAM tum NK NK/ILC Tcm Tcn
1529 1276 213 1587 866 61871 42 422 197 747
Tex Tfh Tha Thn Tp Treg
1969 921 310 3513 358 1071
neighborhoods
cs <- split(colnames(sce), sce$sid)
names(id) <- id <- names(cs)
r <- 0.02; k <- 101
# get radial neighborhood
is <- lapply(id, \(sid) {
sub <- sce[, cs[[sid]]]
xy <- "Center(X|Y)_global_mm"
xy <- grep(xy, names(colData(sub)))
xy <- as.matrix(colData(sub)[xy])
nn <- nn2(xy, searchtype="radius", r=r, k=k)
is <- nn$nn.idx[, -1]; is[is == 0] <- NA; is
})Code
n <- unlist(lapply(is, \(.) rowSums(!is.na(.))))
l <- sprintf("#cells in %sum radial neighborhood", r*1e3)
hist(n, breaks=seq(-.5, max(n)+.5), xlab="", ylab="", main=l)
abline(v=quantile(n, c(.05, .95)), col="red", lw=2, lty=2)
abline(v=m <- mean(n), col="red", lw=2); m[1] 13.55114
composition
# quantify each neighborhood's
# subpopulation composition
fq <- lapply(id, \(.) {
ks <- (se <- sce[, cs[[.]]])$lv2
ks[is.na(ks)] <- "tum"
mx <- matrix(ks[is[[.]]], nrow=ncol(se))
mx[is.na(mx)] <- ""
names(kt) <- kt <- unique(ks)
ns <- sapply(kt, \(k) rowSums(mx == k))
data.frame(
prop.table(ns, 1),
row.names=colnames(se))
}) |> bind_rows(); fq[is.na(fq)] <- 0clustering
set.seed(251210)
km <- kmeans(fq, centers=k <- 10)$cluster
ns <- sprintf("%02d", seq_len(k))
ns <- factor(km, labels=paste0("N", ns))visuals
sce$ctx <- ns[match(colnames(sce), names(ns))]
t(sapply(cs, \(.) table(ns[.]))); table(ns) N01 N02 N03 N04 N05 N06 N07 N08 N09 N10
011LNd 0 2545 0 0 1 0 1306 221 537 56
012LNd 0 4491 0 0 0 0 456 1 45 2
020LNm 13 956 0 0 2 0 929 72 319 4
031LNd 1 4214 0 0 1 0 2227 27 338 8
041SOd 0 2819 0 0 0 0 612 12 120 3
042SOd 0 3000 0 0 0 0 866 29 191 3
051LNn 88 1029 1 0 14 0 2694 314 2612 191
052LNn 147 460 4 0 12 0 1402 626 2486 674
061LNd 1 1437 0 0 1 0 1623 109 740 39
062LNd 0 2344 0 0 0 0 1110 14 198 1
071NPd 752 7 2314 213 12 0 15 117 62 60
072NPd 2287 397 70 0 19 0 1047 578 1502 41
080LNn 0 1892 0 0 2 0 2144 167 697 167
091EYn 0 1413 0 0 0 0 440 30 94 1
092EYn 0 222 0 0 4 0 324 227 338 86
100LNd 1 1122 0 0 4 0 1124 181 676 77
111LNn 114 1049 0 0 37 1 1204 814 906 195
112LNn 20 1110 0 0 110 0 1089 787 532 95
121LNm 5 306 0 0 2525 3435 589 486 282 330
122LNm 12 1725 0 0 765 1214 1942 889 680 861
ns
N01 N02 N03 N04 N05 N06 N07 N08 N09 N10
3441 32538 2389 213 3509 4650 23143 5701 13355 2894